Average Density of States in Disordered Graphene systems 
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In this paper, the average density of states (ADOS) with a binary alloy disorder in disordered 
graphene systems are calculated based on the recursion method. We observe an obvious resonant 
peak caused by interactions with surrounding impurities and an anti-resonance dip in ADOS curves 
near the Dirac point. We also find that the resonance energy (E r ) and the dip position (sdi P ) are 
sensitive to the concentration of disorders (x) and their on-site potentials (v). An linear relation, 
Sdip—xv^ not only holds when the impurity concentration is low but this relation can be further 
extended to high impurity concentration regime with certain constraints. We also calculate the 
ADOS with a finite density of vacancies and compare our results with the previous theoretical 
results. 

PACS numbers: 81.05.Uw, 71.55.-i, 71.23.-k 



I. INTRODUCTION 

Graphene is a two-dimensional (2D) material with 
a single atomic layer of graphite. The material has 
been fabricated firstly by rubbing graphite layers against 
an oxidized silicon surface recently. 1 Due to the linear 
dispersion relation of its electronic spectrum near the 
Dirac point, the electron transport behavior of graphene 
at low energy range is essentially determined by the 
massless relativistic Dirac's equations. Many interest- 
ing properties of graphene have been studied experimen- 
tally and theoretically by many research groups, includ- 
ing the unusual quantum Hall effect quantum min- 
imal conduct ivity^SaLSaS i 10 i 1 1 ferromagnetism, 12,13 and 
superconductivity ! 14 i 15 i 16 

The disorder in graphenes can significantly im- 
pact their electronic properties and has been stud- 
ied extensively^^^^^^^^^ It is believed 
that the interplay between disorders and electron- 
electron interactions determine the low energy behavior 
of the electron in the graphene system. Due to disor- 
ders, the average density states (ADOS) increases at the 
Dirac points. 5 The ADOS is actually an important pa- 
rameter to describe electronic structures, especially in a 
disordered system. The minimum conductivity (cr min ), 
for instance, can be obtained by calculating the diffusion 
constant (D) and ADOS (p) at the Fermi level through 
the Einstein relation (<r m i n = pD). To date, various 
methods have been proposed to calculate ADOS and the 
local density of states (LDOS) in various types of dis- 
ordered graphenes, such as Anderson disorder, 28 short- 
range potential disorder ^ 11 ! 29 ! 30 ! 31 ! 32 , long-range poten- 
tial disorder 11 , and vacancy^* 3 - 3 - However, these calcu- 
lation results are unreliable due to the limitations of ap- 
proximation, or these calculations are only suitable to 
address electronic structure with a low impurity concen- 



tration. It is, therefore, important to obtain the accurate 
electronic structures of disordered graphenes. 

Recently, the phenomenon of spectrum rearrangement 
has been studied in a binary alloy disorder in disordered 
graphene. Due to the limitation of the coherence poten- 
tial approximation (CPA)^ 3 ^ the corresponding results 
are appropriate to address the electronic structures of 
graphene with an extreme low impurity concentration. In 
this paper, we calculate the ADOS with the similar sys- 
tem by using the recursion method. Our numerical sim- 
ulations provide the accurate ADOS for different impu- 
rity concentrations. Moreover, our simulation can also be 
generalized to study other types of disordered graphene. 
Though the main features of ADOS with a binary alloy 
disorder in graphene can be described using the resonant 
and anti-resonant states caused by scattering of a impu- 
rity as reported in Ref. 20. Interestingly, we observe 
that the impurity concentration (x) and the on-site po- 
tential (v) have significant impact on the main features 
of ADOS. The impurity concentration shifts the position 
of resonance energy (E r ). An linear relation for the anti- 
resonance dip shift (edip = xv) with a relative low im- 
purity concentration can be extended to a high impurity 
concentration case under certain conditions. Moreover, 
the ADOS with a finite concentration of vacancies is also 
calculated using our proposed approach and compared 
with the other previous theoretical results. 

This paper is organized as follows. The tight-binding 
model of graphene with a binary alloy disorder is given 
in Sec. II A. The recursion method is introduced and 
its accuracy and applications are explored in Sec. II B. 
The ADOS of graphene with different impurity concen- 
trations and on-site potentials are calculated and the de- 
tailed discussions are given in Sec. III. Finally, we con- 
clude our contributions in Sec. IV. 



2 



II. COMPUTATIONAL MODEL AND METHOD 
A. Model 

Fig. 1 shows the hexagonal lattice structure of a 
graphene and each unit cell has two inequivalent atoms 
labeled A and B, respectively. If we consider the contri- 
bution from the 7r bond (one 7r electron per atom) and 
the nearest interactions in the graphene, the Hamiltonian 
based on the Wannier representation can be expressed as 



H = Y / e i \i){i\-U j Y / \i)(j\, 



(1) 



where, i and j denote the neighboring sites on the lattice, 
Si is the on-site energy, and Uj is the nearest hopping en- 
ergy (Uj=t and its value is close to 2.7 eV in graphenes). 
For simplicity, here t is scaled to be 1. 




FIG. 1: (Color online) A honeycomb lattice of a graphene. 
Here, an unit cell is outlined by the red dashed lines. Two 
inequivalent atoms in the unit cell are labeled with A (the 
green sphere) and B (the blue one), respectively. 

With a given impurity concentration (x), the on-site 
energy Si for a binary alloy disordered system equals v 
with probability x, or zero otherwise. This model, which 
is attributed to Lifshiz, 35 , features the absolutely ran- 
dom distribution of impurities in the space domain. In 
our following calculations, we emphasize on studying the 
properties of the ADOS in two cases. One is on-site po- 
tential (v) varies with given disorder concentrations (x). 
In the other case, the concentration (x) changes with 
given on-site potentials (v). 



B. The Recursion Method 

The Lanczos method 36 (so called the Haydock-Heine- 
Kelly recursion metho d 37 ! 38 ) is a commonly used ap- 
proaches to calculate the ADOS in disorder systems. The 
essential idea of the recursion method is that the Hamil- 
tonian matrix is expressed using a tridiagonal represen- 
tation iteratively. After selecting a localized seed-state 
(|/o))j this recursion method generates a hierarchy of 



states din)) based on a defined orthogonal basis recur- 
sively. 
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where, n = 0, 1, 2, • • • , and the recursive coefficients are 
given by 

a n = {f f f H \ f ;\ b n = ., {fn \ f f n) A bo = 0, !/_!> = 0).(3) 

\Jn\JnJ \Jn-l\jn-l) 

In the orthogonal basis, the Hamilton matrix becomes 



H = 
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(4) 



The diagonal elements in a Green's function matrix 
for a seed state can be derived from Eq. (4) based on the 
continuous-fraction method. 
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and LDOS is defined as 



piocai(E) = lim [ G 00 (E + ie) 

s^0+ 7T 



(6) 



We can then calculate ADOS easily using the following 
equation 



Paver (E) = ^ Plocal(E) , 
M 



(7) 



where M is the number of the samples. To terminate 
continuous fractions, Eq. (5) can be rewritten as 
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In above equation, the terminating term can be written 



t(E) = -{{E - aoo) - [(E - aoo) 2 - 46 2 J}. (9) 

The asymptotic value of the continuous-fraction coef- 
ficient pairs (aoo, b^) in above equation can be obtained 
when n gets large. 

This recursion method has been adopted to investigate 
various disorder systems i 39 i 4Q i 41 In fact, an infinite sys- 
tem can be approximated with periodic boundary condi- 
tion using this method and the numerical error of DOS 
is easy to estimate. 40 If the system size (L) and its recur- 
sion step (TV) are large enough and corresponding posi- 
tive broadening width (e) is reasonably small, the LDOS 
results should be accurate. It is easy to get a reliable 
ADOS of a disorder system by averaging over a large 
number of samples (M). For a perfect graphene, two re- 
lations, e ~ 1/L and p ~ 1/L, hold approximately at the 
Dirac point as shown in Figs. 2 (a) and (6), respectively. 
The results suggest that the DOS in a large finite system 
is close to that in an infinite system. Of course, to ob- 
tain a reliable ADOS result, a large recursion step N is 
needed, which should be larger than L, i.e. 2L [Fig. 2(c)]. 
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FIG. 2: (Color online) (a) £ as a function of L (the blue solid 
line) and the corresponding fitting line shown with red dashed 
line for a perfect graphene. (b) DOS (p) as a function of L 
(the blue solid line) and corresponding fitting line with the 
red dashed line, (c) DOS (p) as a function of N (the blue 
solid line) at the Dirac point (E — 0). Here, L — 800. 

The accuracy of the recursion method is examined by 
comparing the calculated DOS of a perfect graphene with 



TABLE I: The calculated density of state in a perfect 
graphene using both the recursion method and the strict in- 
tegral method within the first BZ. Here, L = 800, N = 2L 
and e = 5x 10~ 3 . 



Energy 


0.0 


0.4 


0.8 


1.2 


1.6 


2.0 


REM a 


0.0037 


0.0781 


0.1990 


0.2570 


0.1952 


0.1702 


INM 6 


0.0037 


0.0781 


0.1998 


0.2529 


0.1938 


0.1695 


Error 


< 0.0001 


< 0.0001 


< 0.002 


< 0.002 


< 0.002 


< 0.002 



a The recursion method 
b The strict integral method 



those obtained based on the strict integral method in 
the first Brillouin zone. The calculated DOS at several 
energy points are list in Table I. Clearly, the results using 
the recursion method are reliable and their relative errors 
are very small. In what follows, we calculate the ADOS 
of a disorder graphene and the parameters, L = 800, 
N = 2L, £ = 5x 10" 3 and M = 1000, are chosen. 



III. RESULTS AND DISCUSSION 

Fig. 3(a) shows several calculated ADOS of graphene 
with the fixed impurity concentration x = 10%, but dif- 
ferent on-site potentials (v). Our results clearly show 
that ADOS remains the similar shape in a high energy 
region. However, the ADOS near the Dirac point changes 
remarkably for different on-site potential (v). Two obvi- 
ous features are observed. (1) An obvious resonance peak 
appears in the ADOS curves when v > 3. Its energy 
position or the resonance energy (E r ) shifts towards the 
Dirac point when the on-site potential (v) increases. This 
peak is attributed to a resonance state. A very similar 
feature was also observed in the LDOS when there is a 
single impurity or a low impurity concentration ] 29 ! 30 (2) 
There is a dip near the Dirac point, which results from an 
antiresonance state. The position of this anti-resonance 
dip (edip) shifts due to the presence of impurities. As 
shown in Fig. 3(b), the linear relation, Sdi P =xv, holds 
when the on-site potential is relatively small, but the dip 
disappears when v becomes larger (i.e.v > 3). These 
results show that the relation (sdi P =xv) is correct when 
v <^ v dip, while this relation does not hold for a system 
with on-site potential v > Vdi P . 

To quantitatively explore the impurity concentration 
how to shift the position of resonance energy (E r ), we 
compare the result for the finite concentration with that 
for a single impurity. The calculated results are shown in 
Fig. 3(c). Here, we adopted the recursion method and the 
effective-mass approximation method (EM A) 42 to deter- 
mined the position of resonance energy (E r ) of LDOS at 
impurity site for a single impurity in the graphene system 
respectively. We also calculate resonance energy (E r ) of 
ADOS for the finite impurity concentration by using the 
recursion method. 

The Green's function in a perfect graphene based on 



4 



the effective-mass approximation is expressed as 
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where S = ?>^/3a 2 /2 is the area of a unit cell in real space. 
a is lattice constant. k c is cutoff wave vector and is set 
to be 2.13/a. 42 For a simple case, the Green's function of 
a single impurity graphene is expressed as 
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LDOS at the impurity site is easy to obtain by 
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The resonant energy (E r ) in this case can be defined 
using the Lifshits equation as 



1 = vReGo(E r ) 



V3v 
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When on-site potential (v) increases as shown in 
Fig. 3(c), the position (E r ) of resonance state shifts to- 
ward to the Dirac point as shown in Fig. 3(a). E r for 
the finite concentration is quite different from that for 
a single impurity. The difference is caused by the mul- 
tiply scattering among impurities. Our simulation im- 
plies that the multiple scattering process is important 
for finite concentration and its contribution neglected in 
some approximation methods is not available. Note that 
the positions of E r for a single impurity determined by 
the recursion method and by the effective approximation 
respectively are still quite difference at a small on-site 
potential and coincide each other at a large on-site po- 
tential. It is because the Green's function obtained by 
Eq. (10) is validate only for low-energy (very close to the 
Dirac point) region. 

Fig. 4(a) presents the calculated ADOS of a graphene 
with different impurity concentrations (x) while the on- 
site potential (v) is set to be 1.5. When the on-site po- 
tential (v) is small, the anti-resonance dip exists, but 
there is no obvious resonance peak. Interestingly, we 
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FIG. 3: (Color online) (a) ADOS (p) of a graphene with im- 
purity concentration x = 10% and various on-site potentials 
(v). (b) The energy shift of the anit-resonance state (edip) 
as a function of the on-site potential (v) (the blue solid line). 
The red dashed line stands for the linear relation, Ed% v — xv. 
(c) The resonance energy (E r ) as a function of the on-site po- 
tential (v). The blue, green, and red lines stand for the ADOS 
for a graphene with impurity concentration x = 10%. ADOS 
is calculated using the recursion method, LDOS of a single 
impurity is computed using the recursion method, and LDOS 
of a single impurity is obtained using EMA, respectively. 



find that the linear relation, Sdi p =xv, still holds even for 
a large impurity concentration, though the multiply scat- 
tering among impurities is strong. Generally, this linear 
relation can be easily understood by using virtual crys- 
tal approximation(VCA) which simply predicts a energy 
level at which two bands coincide is shifted by Sdi^—xv. 
The prediction of VCA is reasonable for a low impurity 
concentration and a small on-site potential, but it fails 
for large values. It is also found that the van Hove sin- 
gularity disappears when the impurity concentration is 
about 40%. This reflects that the phenomenon of the 
spectrum rearrangement happens in the large impurity 
concentration. 

In order to get better understanding of the linear 
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FIG. 4: (Color online) (a) ADOS (p) of a graphene with var- 
ious (number of disorders) impurity concentrations (x). Here, 
on-site potential is set to be 1.5 (v = 1.5). (b) The energy 
shift (sdip) as a function of impurity concentrations (x) (the 
blue solid line), for on-site potential v = 1.5. The red dashed 
line shows the linear relation: Sdiv = xv. 



relation, Sdip—xv^ here, we present a qualitative dis- 
cussion based on the coherent potential approximation 
(CPA). The Self-energy £ can be obtained by neglecting 
multiple scattering process among the impurities as 



xv 



1-vGq(E-H) 



(14) 



The energy shift Sdi v can be determined by substitut- 
ing E — £ =in (k > is real). Thus we have, 
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hi is calculated by the following equation, 
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If <C 1, Eq. (15) can be simplified to edi P =xv. For 
a system with a given impurity concentration e.g. 10%, 



FIG. 5: (Color online) (a) Parameter kv as a function of 
on-site potential (v), where impurity concentrations is x = 
10%. (b) The energy shift (sdip) as the function of on-site 
potential (v). The blue and red solid lines stand for the results 
calculated using the recursion method and CPA, respectively. 
In the calculations, we choose x = 10%. The green dashed 
line indicates the linear relation, Sdip = xv. 



Fig. 6(a) clearly shows that kv is a small value for v less 
than 2.5. If kv <C 1, the relation, edi P =xv, is correct as 
shown in Fig. 3(b). However, the calculated curve of Sdip 
using CPA does not match the linear relation ( K e d i P =xv) 
when v is larger than 2.5. The possible reason is that mul- 
tiple scattering among the impurity cannot be neglected 
in the CPA calculations when v > 2.5. With given the 
impurity concentration 10%, based on the condition of 
\ReE\ = |ira£|, CPA calculation estimates roughly the 
threshold value Vd% v ~ 2.9 of on-site potential (vdi p ) to 
guarantee the linear relation. This result is somewhat 
underestimated by looking at Fig. 3(a), for instance, the 
value of Vdi P should be larger than 3. For the system 
with a given on-site potential (v), Fig. 6(a) shows that 
the results obtained using the CPA method are different 
from those obtained using the recursion method. In the 
case of v = 1.5, the inequality kv <C 1 is satisfied only for 
x < 0.2. With small on-site potential, we find that the 
relation Sdi P =xv still holds even when impurity concen- 
tration is close to 50%. This observation indicates that 
the CPA is not applicable to obtain the reliable results of 
graphene when impurity number is large because which 
neglects the scattering among impurities. 

The threshold impurity concentration (xdip) for v = 
1.5 is about 37% when the liner relationship is correct, 
as show in Fig. 6(b). When the on-site potential (v) 
becomes large enough, Xd% v is extremely small. This re- 
sult is consistent with our calculated results shown in 
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FIG. 6: (Color online) (a) Parameter kv as a function of impu- 
rity concentration (x) when v — 1.5 or v = 5. (b) The energy 
shift (sdip) as a function of impurity concentration (x). The 
blue and dark cyan solid lines stand for the case when v = 1.5 
and v = 5, respectively. The red solid line indicates k as a 
function of impurity concentration (x) for v = 1.5 case. The 
green dashed line shows the linear relation of Sd% v — xv. 



Fig. 5. There is a visible dip only for extreme small im- 
purity concentrations. Clearly, Xdi p is sensitive to the 
on-site potentials. When the inequality raCl does not 
hold anymore, x^ip cannot be obtained through the sim- 
ple CPA calculations i2£ Previous theoretical studie o 3Q i 34 
have predicted that the CPA with single-site scattering 
may fail to predict the dip. Our results obtained by the 
recursion method prove that the CPA is valid only if 
(kv <C 1). The linear relation, Sdi P =xv, can be further 
extended to system with high impurity concentration and 
small on-site potential though the CPA fails to work. 

For large on-site potentials as shown in Fig. 7, there 
are obvious resonance peaks. The height of the peak in- 
creases significantly as the number of disorders increases. 
However, the position of resonance peak does not shift 
obviously at high impurity concentration. These results 
suggest that E r is determined by the on-site potential 
(v) and impurity concentration (x). However, when the 
on-site potential is large, E r does not depend obviously 
on the impurity concentration even x has large value. It 
is also found that the van Hove singularity disappears 
when the impurity concentration is about 20%. 

The disordered graphene with a finite density of vacan- 
cies can be modeled by setting the on-site potential v to a 
very large value i 32 i 42 Fig. 8 shows the calculated ADOS of 
graphene systems with vacancies (v = 1000) based on the 
recursion method. There is a clear sharp peak near the 
Dirac point, which is in line with the numerical results ob- 



tained using the stochastic recusive method^This sharp 
peak can be fit very well by the Lorentz distribution. 
However, the precious theoretical calculation, based on 
the full Born approximation (FBA), 5 predicted the value 
of ADOS should be exactly zero. Its error may occur be- 
cause the multiple scattering is ignored in the FBA cal- 
culation. Meanwhile the calculation based on the CPA 4 
and the full self-consistent Born approximation (FSB A) 5 
also fails to produce observable resonance peaks. These 
results indicate that there are some limitations in the 
CPA, FBA and FSBA methods near the Dirac point. 
The correctness or accuracy of these methods need to be 
further examined. 




FIG. 8: (Color online) The calculated ADOS (p) of graphene 
with different concentrations of vacancies. Here, we chose 
v 1000. 
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IV. CONCLUSION 

In this paper, these ADOS of a binary alloy disorder 
in disordered graphene are calculated using the recursion 
method. The applicability and accuracy of the recursion 
method are addressed. It is found that the shape of the 
resonance peak and the position of the anti-resonance 
dip are sensitive to impurity concentration (x) and on- 
site potential (v). The linear relation, Sdi P =xv, can be 
derived based on the CPA when ctCI. This relation is 
able to explain the position shift of anti-resonance dips at 
low or high impurity concentration for small on-site po- 
tential case. For large u, either the CPA or Eq. (14) 
fails except for extremely low impurity concentration 
(kv <C 1). The main reason for the failure of the CPA 
is that the scattering among impurities is neglected. By 



setting v to be a huge value, the model can be used to 
simulates finite concentration of vacancies in a graphene. 
The resonance peak of ADOS at the dirac point is found. 
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